/*
 *  viscoelastic.h
 *  EPPI-0.0
 *
 *  Created by Sergio Conde Martin on 1/04/11.
 *  Copyright 2007 __MyCompanyName__. All rights reserved.
 *
 */
 
#ifndef _viscoelastic_h_
#define _viscoelastic_h_

#include <iostream>
#include <sstream>
#include "MATERIALS/material.h"

using namespace std;

class Viscoelastic : public Material
{
public:
		Material *matHyper;
private:
		int nTau;	 				 // Number of internal variable
		vector < double > vTau;		 // Vector of relaxation (or retardation) time, , as many as number of internal variables.
		vector < double > vBeta;	 // Vector of strain-energy factor, as many as number of internal variables.
		VectorXd vMu;
		VectorXd vEta;

		bool integrated;				// Pre-integrated or integrated? Integrated == 0 >> I. Variables to be integrated. Integrated == 1 >> I. Variables pre-integrated
        bool rheological;            // Is Viscoelasticity rheological?
        double sumbeta;
        VectorXd Idv;
        MatrixXd Id4m;

        Tensor2Sym mSTens, mETens, GamTens;
        Tensor4LSym CtrheTens;
	
        //*** Consistent Integration ***

        double Wrhe;
        VectorXd Srhe, E, Ci, tildeGam;
        MatrixXd Ctrhe, D, PMat, PtMat, tildePMat, tildeGamoCi, CiotildeGam;
	
public:
        Viscoelastic();
        Viscoelastic(ifstream &, const string &, const string &,  Material*);
        ~Viscoelastic(){}
	
        /** Virtual functions from base */

        void calcW(const VectorXd &, const VectorXd & );	// Computes the viscoelastic free energy function
        void calcS( const VectorXd &, const VectorXd & );	// Computes the viscoelastic Second Piola-Kirchoff Stress Tensor
        void calcQv(const int &, const VectorXd &, const VectorXd & );	// Computes the Non-Equilibrium Stress
        void calcCt(const VectorXd &, const VectorXd & );             // Computes the constitutive material tensor


        int getNIntVars(){return nTau;};
        bool getMode(){return integrated;};
	
        VectorXd getTau();
        VectorXd getBeta();
        VectorXd getEta(){return vEta;};
        VectorXd getVMu(){return vMu;};

private:
        /**Local Function */
        void calcWrhe( const int &, const VectorXd &);
        void calcSrhe( const int &, const VectorXd &);
        void calcCtrhe( const int &, const VectorXd &);
    
        void calcDcPsiGamma (const int &, const VectorXd &, const VectorXd &);
};

#endif
